clear
%oad ../parameter_sandbox_ces/eq_mit_mass.mat 
load ../../prod/calibration_v3.5/ces/eq_mit_mass.mat
eq_mit_ces = eq_mit;

%load ../parameter_sandbox_wage_rule/eq_mit_mass.mat 
load ../../prod/calibration_v3.5/kimball/eq_mit_mass.mat

%%
Ns = 1 - [0 6.69652130]/100;
 prods = repmat(exp(glob.sf(:,2)), 1, options.T);

%% pick scaling coefficients
addpath ../
adj = sum(eq_mit_ces.mu,1);

scx_ces = pick_scale(adj, Ns);
%scx_ces = [0 -diff(Ns-1)]/scx;

adj = sum(eq_mit.mu,1);

scx = pick_scale(adj, Ns);
%scx = scx_ces;%.0645/scx;

%% compute cost-weighted markups
ngdp = sum(eq_mit.mu.*eq_mit.p.*eq_mit.y);
labor_bill = eq_mit.W'.*sum(eq_mit.l.*eq_mit.mu);
labor_share = (labor_bill./ngdp);

ngdp = sum(eq_mit_ces.mu.*eq_mit_ces.p.*eq_mit_ces.y);
labor_bill = eq_mit_ces.W'.*sum(eq_mit_ces.l.*eq_mit_ces.mu);
labor_share_ces = (labor_bill./ngdp);


%%
set(0,'defaultaxesfontname','cambria math') % beautify the axes a bit
set(0,'defaultTextFontName', 'cambria math')


close all

mu = (log_linear_scale(1./labor_share, scx));
mu_ces = (log_linear_scale(1./labor_share_ces, scx_ces));

Ys    = repmat(eq_mit.C',glob.Nsf,1);
Z     = sum(eq_mit.y./Ys./prods.*eq_mit.mu);
A   = log_linear_scale(Z.^(-1), scx);

Ys    = repmat(eq_mit_ces.C',glob.Nsf,1);
Z     = sum(eq_mit_ces.y./Ys./prods.*eq_mit_ces.mu);
A_ces   = log_linear_scale(Z.^(-1), scx_ces);


%%
bl = three_equation_model(A, mu, .5);

ces_mu = three_equation_model(A, mu_ces, .5);
ces_tfp = three_equation_model(A_ces, mu, .5);

ces_bl = three_equation_model(A_ces,mu_ces, .5);

%% plot
close all
subplot(2,1,1)
plot(bl.Y, 'LineWidth', 4)
hold on
plot(ces_tfp.Y,'--', 'LineWidth', 4)
plot(ces_bl.Y, '-.','LineWidth', 4)
legend('Baseline', 'CES + Markup variation', 'CES', 'Location', 'SouthEast')
xlim([1 T])
title('Y')

subplot(2,1,2)
plot(bl.L, 'LineWidth', 4)
hold on
plot(ces_tfp.L,'--', 'LineWidth', 4)
plot(ces_bl.L,'-.', 'LineWidth', 4)
xlim([1 T])
title('L')

% subplot(2,2,3)
% plot(bl.W, 'LineWidth', 4)
% hold on
% plot(ces_tfp.W,'--', 'LineWidth', 4)
% plot(ces_bl.W,'-.', 'LineWidth', 4)
% xlim([1 T])
% title('W')
% 
% subplot(2,2,4)
% plot(mu_ces, 'LineWidth',4)
% hold on
% plot(mu,'--','LineWidth',4)
% xlim([1 T])

set(gcf,'units','points','position',[10,10,500,500])
set(findall(gcf,'-property','FontSize'),'FontSize',16)

print('-dpng', 'kimball_ces_diff.png')